install.packages("devtools") 
devtools::install_github("jonathandroth/pretrends")

library(readr)
library(tidyr)
library(pretrends)

##Firm-level COD Emission

#Below College
#Input beta
data <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_belowcl_beta.csv")
data_long <- pivot_longer(data, cols = everything(),values_to = "value")
beta <- data_long$value

#Input sigma
sigma <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_belowcl_sigma.csv")

#t
tVec<-c(-5,-4,-3,-2,0,1,2,3,4)
referencePeriod <- -1 #This is the omitted period in the regression

data.frame(t = tVec, beta = beta)

#Compute slope that gives us 90% power
slope80 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.8,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope90 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.9,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope80

slope90

#Belowhs
#Input beta
data <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_belowhs_beta.csv")
data_long <- pivot_longer(data, cols = everything(),values_to = "value")
beta <- data_long$value

#Input sigma
sigma <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_belowhs_sigma.csv")

#t
tVec<-c(-5,-4,-3,-2,0,1,2,3,4)
referencePeriod <- -1 #This is the omitted period in the regression

data.frame(t = tVec, beta = beta)

#Compute slope that gives us 90% power
slope80 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.8,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope90 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.9,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope80

slope90

#Rural
#Input beta
data <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_rural_beta.csv")
data_long <- pivot_longer(data, cols = everything(),values_to = "value")
beta <- data_long$value

#Input sigma
sigma <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_rural_sigma.csv")

#t
tVec<-c(-5,-4,-3,-2,0,1,2,3,4)
referencePeriod <- -1 #This is the omitted period in the regression

data.frame(t = tVec, beta = beta)

#Compute slope that gives us 90% power
slope80 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.8,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope90 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.9,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope80

slope90

#Skilled Occupation
#Input beta
data <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_occupation_beta.csv")
data_long <- pivot_longer(data, cols = everything(),values_to = "value")
beta <- data_long$value

#Input sigma
sigma <- read_csv("G:\\Dropbox\\Environmental Injustice\\Empirical\\results\\Table S14_Firm_occupation_sigma.csv")

#t
tVec<-c(-5,-4,-3,-2,0,1,2,3,4)
referencePeriod <- -1 #This is the omitted period in the regression

data.frame(t = tVec, beta = beta)

#Compute slope that gives us 90% power
slope80 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.8,
                  tVec = tVec,
                  referencePeriod = referencePeriod)

slope90 <-
  slope_for_power(sigma = sigma,
                  targetPower = 0.9,
                  tVec = tVec,
                  referencePeriod = referencePeriod)


slope80

slope90